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OBJECTIVES 

The purpose of this study is to establish the technical ground for modeling 
the physics of laser powered pulse detonation phenomenon. Laser powered 
propulsion systems involve complex fluid dynamics, thermodynamics and 
radiative transfer processes. Successful predictions of the performance of laser 
powered launch vehicle concepts depend on the sophisticate models that reflects 
the underlying flow physics including the laser ray tracing the focusing, inverse 
Bremsstrahlung (IB) effects, finite-rate air chemistry, thermal non-equilibrium, 
plasma radiation and detonation wave propagation, etc. The proposed work will 
extend the base-line numerical model to an efficient design analysis tool. The 
proposed model is suitable for 3-D analysis using parallel computing methods. 
The technical objectives of the present project are: 

1 . Extend the laser propulsion analysis model to 3-D applications. 

2. Develop the laser ray-tracing model and laser-focusing model suitable for 
general 3-D geometry conditions. 

3. Implement the efficient treatment for modeling the plasma radiation model 
such that it can be ready for 3-D applications. 

4. Investigate the laser ignition process and its effects on the launch vehicle 
performance. 

5. Develop the ray-tracing and radiation model to be suitable for parallel 
computing environment. 

6. Investigate the air-breathing and rocket modes modeling effects. 
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INTRODUCTION 


Currently, NASA’s aim of operating low cost launch and space vehicles 
requires the research and development of advanced propulsion technologies and 
concepts. One plausible advanced concept is the utilization of off-board pulsed 
laser power source to propel small payload (e.g. 100kg) into earth orbit. The 
merit of the laser-propelled vehicles is in its high efficiency (do not need to carry 
fuel) and high specific impulse. Previous SDIO research led to the invention of 
the one of the laser powered launch vehicle concept - the Laser Lightcraft 
concept, currently being tested at the High Energy Laser Test System Facility, 
White Sands Missile Range, New Mexico. Although the spin-stabilized small 
scale Lightcraft model (invented by Myrabo) has been flown successfully up to an 
altitude of 30 meters using a 10 kW pulsed-laser at 10 Hz, many technical issues 
need to be addressed before an optimized design of the vehicle and its operation 
can be achieved. 

The purpose of this study is to establish the technical ground for modeling 
the physics of laser powered pulse detonation phenomenon. The principle of the 
laser power propulsion is that when high-powered laser is focused at a small 
area near the surface of a thruster, the intense energy causes the electrical 
breakdown of the working fluid (e.g. air) and forming high speed plasma (known 
as the inverse Bremsstrahlung, IB, effect). The intense heat and high pressure 
created in the plasma consequently causes the surrounding to heat up and 
expand until the thrust producing shock waves are formed. This complex 
process of gas ionization, increase in radiation absorption and the forming of 
plasma and shock waves will be investigated in the development of the present 
numerical model. In the first phase of this study, laser light focusing, radiative 
absorption and shock wave propagation over the entire pulsed cycle are 
modeled. The model geometry and test conditions of known benchmark 
experiments such as those in Myrabo’s experiment will be employed in the 
numerical model validation simulations. The calculated performance data (e.g. 
coupling coefficients) will be compared to the test data. Plans for the numerical 
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modeling of the detailed IB effect will also be described in the proposed 
investigation. The final goal will be the design analysis of the full-scale laser 
propelled flight vehicle using the present numerical model. 

In this study, Engineering Sciences, Inc. has developed a laser powered 
launch vehicle performance analysis tool based on its in-house flow and radiation 
codes. UNIC-UNS unstructured-grid flow code and GRADP-UNS unstructured- 
grid radiation code are two advanced numerical models. Many complex 
engineering design problems related to fluid dynamics and radiative heat transfer 
have been solved using these two codes. However, high-temperature 
thermodynamics and plasma dynamics models still need more developments 
and benchmark data validations before these two codes can be applied to the 
design problems of laser powered launch vehicles. The development work shall 
include transient shock capturing algorithm using unstructured-grid method with 
dynamic local refinement and coarsening. High temperature thermodynamics 
and plasma gas dynamics physics shall be modeled and validated. Non- 
equilibrium radiation model with the effects of gas breakdown and laser energy 
absorption also needs to be addressed and modeled. These advanced 
thermodynamics and radiation models will serve as the fundamental building 
blocks for the development and utilization of the proposed laser powered launch 
vehicle performance prediction and analysis tool. 

As a continuing effort of improving the capability of the laser propulsion 
modeling analysis tools, general formulations and modeling specifics outlined in 
the previous section are planned and implemented in the present project. In the 
first part of this project, a three-temperature model is implemented in the baseline 
CFD code [1] and tested for the model-A laser lightcraft configuration at 400 J 
power level. In this study, it is found that the previously employed model for 
evaluating the Landau-Teller relaxation time scale is not appropriate for three- 
temperature calculation. A more general model [2] (Peter A. Gnoffo, et al., 1989) 
is therefore implemented in the present study. Since the three temperatures 
(namely the heavy-particle temperature, vibrational temperature and electron 
temperature) are strongly coupled in the transient process, an implicit solution 
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method is employed to provide stable solutions for this system. The present 
three-temperature model is described in the following section. 

For testing the model performance, the three-temperature model is 
employed for the thrust performance computations of the model-A laser lightcraft 
configuration at 200 J, 300 J, 400 J, 600 J and 800 J power levels. The plasma 
resonant energy conversion factor is fixed at 0.4. Another approach using an 
idealized model described by C. S. Liu, et al. [3], has been investigated. The 
goal of this investigation is to replace the assumed plasma resonant energy 
factor with a better approach. In Liu’s model, the absorption coefficient of a 
plasma sheet beyond resonance frequency can be estimated as a function of 
incoming ray angles and a critical length scale. Preliminary results of coupling 
coefficient calculation using this model have shown under prediction in vehicle 
thrust. This indicates that the conversion factor predicted is on the low side. The 
coupling coefficient predicted using Liu’s model is too low compared to the test 
data (about 19.5% for the 400 J case). It is concluded that Liu’s model, which is 
derived using many simplified assumptions, is not general enough for all laser 
incident angles. The low absorption level predicted using this model is not 
acceptable to be incorporated in the current method. Therefore, it is decided to 
stay with the current fixed 40% absorption rate for plasma resonance. 

Extension of the computational model to 3-D operation conditions using 
parallel computing methods is also completed in the present effort. The 3-D 
model involves more complicated general ray tracing and boundary reflection 
conditions for the radiation model. Parallel algorithms for the general ray tracing 
and radiation solution are implemented in the present analysis tool. Some 
observations and recommendations for the air-breathing and rocket mode 
applications are discussed. 
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COMPUTATIONAL FLUID DYNAMICS MODEL 


Governing Equations 

For The Continuity, Navier-Stokes and Energy (Total Enthalpy) Equations, 
can be written in a Cartesian tensor form: 
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where p is the fluid density, u, is the f h Cartesian component of the velocity, p is 
the static pressure, p is the fluid viscosity, P r is the Prandtl number, H is the gas 
total enthalpy and V stands for the sum of velocity squared. In Eq. (4), k b , n g , T e , 
X e , Q r , Qv and Q ec are the Boltzmann’s constant, electron number density, 
electron temperature, electron thermal conductivity, radiative heat source from 
laser absorption and radiative transfer, vibrational-translation energy transfer 
source term, and the energy transfer due to eletron/particle elastic collisions, 
respectively. The shear stress can be expressed as: 
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The species conservation equation is expressed as: 
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where Y\ is the f h species mass fraction, D is the mass diffusivity, c Y is the 
turbulent Schmidt number, and (b i is the chemical reaction rate for species /' 
respectively. 

Vibrational Energy Equation 

For high temperature gas flows, thermal non-equilibrium state may be 
important. In Landau and Teller’s derivation, a master equation is employed to 
describe the evolution of the population of quantum level N|. This master 
equation is written as: 


dN : 


, L = n 2 K j~ n j~ n L N i ; 1 = O’ 1 ’ 2 

dt ;=0 j=0 


Results from the quantum mechanical solution of the harmonic oscillator 
are used to relate the various quantum transition rates to one another, and then 
the master equation may be summed over all quantum states to arrive at the 
Landau-Teller equation: 
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where p, e v , e v eq and t LT represent the gas density, vibrational energy, effective 
(equilibrium) vibrational energy and the vibrational-translational relaxation time 
scale respectively. An empirical expression (to be discussed in the next section) 
is used to model the Landau-Teller relaxation time scale. 
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As discussed by Gnoffo (1989), the vibrational-translation energy 
relaxation time scale can be evaluated using the following equation 

s=nu>l 

s=mol. 


and 
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where subscript s represents the participating species (only diatomic species are 
involved here), a, n, m and // denote the effective cross section, species number 
density, species unit mass and reduced molecular weight, respectively. The 
correlation constant, A s , has different values for each participating species. They 
are given in the table below. 
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BENCHMARK TEST CASES STUDY 


400 J 30-usec Pulse Width Test Case 

The test conditions of this case are listed below. The laser ignition 
procedure of the baseline model is used. A pie-area of spark ignition region is 
assigned initially. It is assumed that 15 percent of the laser energy is absorbed 
and converted into thermal energy during this spark ignition time period. When 
laser absorption efficiency reaches 1 5 percent, the spark ignition process will be 
terminated and normal laser absorption mode would remain. After this point, all 
absorbed laser energy is used as the energy source for the electron energy 
equation. Subsequent energy transfer to heavy particle and vibrational mode is 
solved by the coupled three-temperature energy equation. When the electron 
number density reaches the threshold value beyond which plasma resonance 
would occur, 30 percent energy conversion rate is assumed for the present 
model. More rigorous model should be employed to model the energy conversion 
for plasma resonance. The flowfield solutions at two time levels are shown in the 
following figures (Figures 1 to 4). The coupling coefficient time histories of two- 
temperature and three-temperature models are also compared in Figure 5. 

Model Conditions: 

• Vehicle Configuration Model-A of Myrabo’s Flight Tests (Closed Air Inlet) 

• Pulsed Laser at 1 0 Hz and 30 ps Pulse Width 

• Laser Energy Level = 400 J 

• Total Number of Mesh Elements = 16,870 

• Assuming Axisymmetric Model 

• Assumed Energy Absorption Scheme 

• Testing the Predictions of Detonation Wave Propagation and Vehicle Thrust 
Calculation 
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Figure 1 . Pressure contours (ATM) at 300 time steps. 
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Figure 3. Pressure contours (ATM) at 1 170 time step. 
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Figure 5. Comparisons of coupling coefficient time histories for two-temperature 

and three-temperature models. 

This three-temperature model, with 30 percent absorption rate assumed 
for the resonant plasma, under-predicts the coupling coefficient as shown in 
Figure 5. The result of the two-temperature model is very close to the measured 
data. This is expected since the assumed energy conversion rate in the plasma 
resonance regime of 30 percent was tuned based on the two-temperature model. 
A higher energy conservation rate will give good predictions for the current three- 
temperature model. A more rigorous model will be employed to eliminate this 
energy conversion assumption. 

Resonant Energy Conversion Factor Study 

Model-A lightcraft thrust coupling coefficients for 400 J, 600 J and 800 J 
power levels have been predicted using the present 3 temperature computational 
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model. The computation starts with a small time step size, 1 .OE-7 sec, and steps 
up to larger time steps within the laser pulse width of 30 micro-sec. Small time 
step at the start is necessary for stable solution due to high concentration of laser 
energy near the focus. A typical case takes about 8000 time steps to produce a 
nearly constant integrated coupling coefficient. The grid size for this model is 
17,184. It takes around 3 days run time to obtain a solution for each case. Figure 
6 shows the coupling coefficient comparisons for 400 J laser power level. 



Figure 6. Comparisons of coupling coefficient time histories for two-temperature 
and three-temperature models for 400 J laser. 

Figures 7 and 8 show the predicted coupling coefficients for 600 J and 800 
J laser power level. 
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COUPLING COEFF. ® COUPLING COEFF. 



7. Predicted coupling coefficient for 600 J laser power. 


Figure 
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The present three-temperature model over-predicts the coupling 
coefficient slightly. Lower coupling coefficient prediction is expected as non- 
equilibrium plasma radiation model is considered in the computation. 


Study of All Power Levels 

Model-A lightcraft thrust coupling coefficients for 200 J, 300 J, 400 J, 600 
J and 800 J power levels have been predicted using the present 3 temperature 
computational model. The computation starts with a small time step size, 1 .OE-7 
sec, and steps up to larger time steps beyond the laser pulse width of 30 micro- 
sec. Small time step at the start is necessary for better time resolution and 
solution stability due to high concentration of laser energy near the focus. A 
typical case takes about 8000 time steps to produce a nearly constant integrated 
coupling coefficient. The grid size for this model is 17,184. It takes around 3 
days run time to obtain a solution for each case using one processor. Figure 9 
shows the coupling coefficient comparisons for the 5 laser power levels. 

In the present spark ignition modeling, all spark energy (which is assumed 
to be 1 5% of the available beam power) is used as the source of the heavy gas 
energy equation. Within the spark time period, the one-temperature model is 
used. Vibrational and electron temperatures are then assigned to be the same 
level as the heavy gas temperature after every time step. This treatment is found 
to create inconsistent energy source budget when different initial time step size is 
used. The three-temperature model seems to produce consistent results for 
different initial time step sizes. This approach will be further investigated to 
confirm that it will be the better alternative than the present spark ignition model. 
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Figure 9. Comparisons of coupling coefficient time histories for 5 laser power 

levels. 


Parallel and 3-D Rav Tracing Radiation Model 

a. Laser Beam Rav T racing Using Parallel Computation 

Laser radiation absorption is one of the most important phenomena in the 
analysis of laser propelled lightcraft and its modeling is based on the geometric 
optics in the present study. In the geometric optics, the laser beam is assumed to 
consist of a finite number of individual laser rays, each ray is then traced until it is 
totally absorbed by a gas or surface element. In the previous analysis, ray tracing 
for the laser beam was implemented on a sequence computer. In the current 
period of the project, a parallel algorithm has been developed to extend the ray 
tracing calculation of laser beam in parallel computation. 
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The ray tracing is a time sequence process and its implementation in a 
sequence computer is straightforward through tracing a ray volume by volume 
until a boundary is reached. However, when multiple computers are used in a ray 
tracing process, its implementation becomes different since only one computer 
(not all computers) is involved in an actual ray tracing calculation while other 
computers are essentially idle. To implement a parallel algorithm, an iterative 
computation procedure is developed in the present computer program. Each 
iteration involves ray tracing in a single subdomain associated with a specific 
computer. In a specific iteration, all information required in the ray tracing is 
passed from the previous iteration. A detailed implementation procedure is 
provided below. 

(1) Decompose the computational domain by the METIS software. Each 
subdomain is then assigned to a specific computer processor. 

(2) Assign the ray tracing flag values to all the computer processor. The value 
of 1 means the ray tracing takes place in the current processor. 

(3) Decide the subdomain and the corresponding processor index number 
from which each ray starts for a ray tracing calculation. The value of the ray 
tracing flag is set to be unity for this processor. 

(4) Start the ray tracing in the processor with flag= 1 while other processors are 
idle. The ray tracing calculation involves tracing a ray volume by volume in 
the designated subdomain until a boundary element is reach. 

(5) If the boundary element is a zonal boundary, the information on the 
intersection point is saved and passed to its adjoint subdomain by calling 
the PVM library. The ray tracing flag becomes unity for the processor 
corresponding to this adjoint subdomain. Then, the calculation enters into 
the next iteration by repeating step (4). If the boundary element is reflective 
surface, the ray tracing continues in the same subdomain after the ray 
reflection. If the boundary element is other surface, the ray tracing ends 
here and the calculation move on for the next ray. 
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In order to check the accuracy and efficiency of the developed parallel 
algorithm, a laser propelled lightcraft problem with the laser power equal to 370 J 
has been examined for three cases with processor number N p equal to 1, 2, and 
6, respectively. Figure 10 demonstrates the corresponding coupling coefficient 
distributions. It is clear that the results from different processor number are very 
close. This confirms that the present parallel implementation is correct. In 
addition, the parallel performance is also investigated and the speedups for 
different N p are found to be close to the ideal speedup. 



Time, ms 

Figure 1 0. Comparison of the coupling coefficient distributions for different 

processor number. 

b. Laser Beam Rav Tracing in 3-D Geometry 

In the previous analysis, the ray tracing calculation was only applied for a 
2-D geometry. In the current period of the project, the computer code has been 
developed to extend the ray tracing calculation for a 3-D geometry. In the 2-D ray 
tracing code, the laser ray tracing conditions are easily provided in the input data 
file unic.inp. However, the ray tracing conditions for a 3-D geometry becomes 
much more involved. Thus a special subroutine has been developed to 
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implement these conditions in the present 3-D study. In this subroutine, the 
geometrical and laser intensity information for the input laser imprint is provided. 
The imprint area is then discretized into many small surfaces. Each of these 
surfaces is associated with a laser ray. It is noted that the descretization of the 
laser imprint area has nothing to do with the computational grid used for the laser 
lightcraft. To connect these two grids, the imprint area grid along with the laser 
ray intensity is mapped on the boundary of the computational grid. With this 
mapping, the starting location and power intensity for each ray on the 
computational grid are known and the ray tracing calculation can begin. Similar to 
a 2-D geometry, the ray tracing in a 3-D geometry involves tracing a ray volume 
by volume until a boundary is reached. 

In order to validate the present 3-D ray tracing code, a 3-D problem has 
been constructed by rotating the existing 2-D grid over 45°. The laser imprint 
area considered is of square shape. Due to the symmetry, only a eighth of the 
square is used. The imprint area is divided into 30x20 grid representing 600 
laser-beam rays. By analyzing all the intersection points of each ray in the 
computational domain, the 3-D ray tracing has been found to be correct. With the 
validation of 3-D ray tracing calculation, a 3-D modeling of laser propelled 
lightcraft becomes possible and it has been tried in a hypothetical 3-D problem. 

Air Breathing and Rocket Mode Effects 

In practical design of laser powered launch vehicles, the performance of 
air supply through intake passages may affect the overall thrust performance of 
the vehicle. At very high altitudes, it will be necessary to switch to rocket mode 
that the propellant carried in the vehicle are injected and detonated for each 
pulsed laser cycle. For numerical analysis for the air breathing or rocket modes, 
the same physical models can be applied. What needs to be included in the 
model is the inclusion of air intake passage control and injection scheme 
arrangement. It is envisioned that the present model can be directly applied to 
these conditions without any problem. 
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CONCLUSIONS 


In the present research project, a comprehensive computational model 
has been developed for the prediction of laser supported propulsion system that 
involved thermal, chemical and radiative non-equilibrium flow physics. The 
model-A laser lightcraft design has been analyzed for a series of laser power 
inputs. The predicted thrust coupling coefficients show good agreement with the 
measured data. Extension of this model to 3-D and parallel computing 
conditions, which involves a general radiative ray-tracing scheme, has been 
completed and tested. More studies are needed to determine the plasma- 
resonance energy conversion factor, which seems to be of major contributions to 
the integrated thrust coupling coefficients. 
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